Instability of current sheets and formation of plasmoid chains 
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Current sheets formed in magnetic reconnection events are found to be unstable to high- 
wavenumber perturbations. The instability is very fast: its maximum growth rate scales as 
S 1/I4 va/Lcs, where Lcs is the length of the sheet, va the Alfven speed and S the Lundquist number. 
As a result, a chain of plasmoids (secondary islands) is formed, whose number scales as S 3 ^ 8 . 

PACS numbers: 52.35.Vd, 52.35.Py, 94.30.cp, 96.60.lv 



Magnetic reconnection is a plasma phenomenon in 
which oppositely directed magnetic field lines are driven 
together, break and rejoin in a topologically different con- 
figuration. It is an essential element in our understanding 
of the solar flares and the magnetotail, 0,0,0] where it is 
directly observed, [J, Q as well as of other astrophysical 
plasmas. On Earth, it plays a crucial role in the dynam- 
ics of magnetically confined plasmas in fusion devices. @ 
There are two standard reconnection models: the 
Sweet-Parker (SP) model, @, § and the Petscheck 
model. @ The latter is very appealing as it predicts fast 
reconnection rates similar to the observed ones. However, 
numerical simulations have consistently failed to repro- 
duce it unless a spatially inhomogcneous anomalous re- 
sistivity is usedfl(| or Hall physics is invoked. 0, [ll[ In 
contrast, the SP reconnection, characterized by long cur- 
rent sheets (~ system size) and slow reconnection rates 
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where 77 is the plasma resistivity, is routinely ob- 
served both in experiments [l2| and in simulations. 

The break up of current sheets and formation of plas- 
moids (secondary islands) appears to be a generic fea- 
ture of reconnecting systems. Plasmoids have been ob- 
served both in solar flares [HI and in the magnetotail. 14 
In numerical simulations, plasmoid formation has been 
reported in many different set ups, from fluid (l5l fl6j to 
fully kinetic. [U 



18| Plasmoids have been popular in the- 
ories of magnetic reconnection and related phenomena: 
e.g., they have been invoked as a plausible mechanism 
for accelerating reconnection, cither by decreasing the 
effective length of the SP current sheet and/or by trig- 
gering anomalous-resistivit y m echanisms associated with 
small-scale plasma effects :|19| a multiple-plasmoid sce- 
nario has been suggested to explain the production of 
energetic electrons during reconnection events: [20j| it has 
been conjectured that periodic ejection of plasmoids in 
star-disk systems could account for the knot-like struc- 
tures observed in stellar iets-plj 

A theoretical understanding of the mechanism whereby 
the plasmoids arc formed has, however, been lacking. 
It has been believed that plasmoid formation is due to 
a standard tearing instability [22j of the current sheet. 
Bulanov et al. 23:1 considered the onset of such an in- 



FIG. 1: (Color online) Schematic drawing of the current sheet 
with in- and outflows and a forming plasmoid chain. 



stability by taking into account a linear outflow along 
the sheet and repeating the standard tearing-mode cal- 
culation. While this gave a qualitative indication of the 
mildly stabilizing role of the outflow, it was necessarily 
a nonrigorous approach because the resistivity could not 
be neglected anywhere inside the current sheet and no 
other small parameter was available. In this Letter, we 
pursue a different strategy by considering a large-aspect- 
ratio sheet and using the inverse of its aspect ratio as 
the small parameter. We show that m ultip le plasmoids 
form and that, unlike in Bulanov et aL.|23l| the instabil- 
ity is very fast, with a rate much larger than the ideal 
(Alfvenic) rate. 

Theoretical estimates for the speed of magnetic recon- 
nection in natural systems arc usually based on the idea 
that a current sheet is formed (Fig. [1]), with the length 
Lcs determined by the system's global properties and 
the width Lcs/v^, where S = vaLcs/v is the Lundquist 
number, Va is the upstream Alfven speed, and r\ is the 
magnetic diffusivity. Let us consider such a current sheet, 
and ask if it can be stable over any significant period of 
time. The SP reconnection time is {Lcs/ v a)vS. If we 
consider times much shorter than this time, we can de- 
termine the structure of the current sheet by seeking a 
stationary solution of the induction equation 



d t B + u ■ VJ3 = B Vu BV ■ u + i]\7 2 B, 



(1) 



where u is the velocity field and B the magnetic 
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field, which we will measure in velocity units. While 
the desired resistive equilibrium is stationary, it is not 
static. It is a consistent feature of current sheets, both 
measured (l2j and simulated, that they support lin- 
ear outflows along themselves. If incompressibility is as- 
sumed and the reconnection is considered in two dimen- 
sions (x, y) with a current sheet along the y axis, then, 
inside the current sheet, u x = — Tqx and u y = Toy. Here 
r = 2va/Lqs, so the outflows are Alfvenic. Obviously, 
this choice of u is a drastic simplification of the real flow 
profile (ignoring, for example, the flow vorticity). It is 
only intended as a minimal model incorporating what we 
believe to be the key feature, namely the linearly increas- 
ing Alfvenic outflow. 

A simple equilibrium solution Bo of Eq. ([1]) exists 
that accommodates this flow pattern: Bq x = 0, and 
Bq v = Bq v (x) satisfies 




FIG. 2: (Color online) The equilibrium magnetic field Bo y 
and the current density joz = d x Bo y (normalized to va and 
va/Scs, respectively). 



d x (xB 0y ) = 0, 



(2) 



where Sqs = (v/^o) 1 ^ 2 is the characteristic width of the 
current sheet. The solution of this equation that vanishes 
and changes sign at x = (the center of the sheet) is 
B 0y = VAf(0) where £ = x/S C s and 
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= ae w 



dze z2 / 2 . 



(3) 



The integration constant a is chosen by matching this so- 
lution with the magnetic field outside the sheet: Bq v = 
ivA, which is a solution of Eq. ([T]) with constant inflows 
u x = T u o, u y = 0. In order to complete our simple 
model of the current sheet, we must choose a suitable 
point x = ixo at which the inside and outside solu- 
tions could be matched. The flow is discontinuous at 
this point: u = (uo,0) for x < — xq, u = (— ToX,Toy) 
for — xq < x < xq and u = (— uo,0) for x > xq, 
where uq — TqXq. We consider this to be an accept- 
able simplification of the real, more complicated, and 
continuous, flow profile. The natural matching point is 
the point where Bo y (x) has its maximum (minimum), 
Boy(±xo) = 0, and where, therefore, the current van- 
ishes. This gives xo = £o<$cSi where £o « 1.31. We then 
require /(±£q) = ±1, so a = £o- The equilibrium mag- 
netic field and current are plotted in Fig. [5] 

We shall show that this equilibrium is subject 
to a very fast linear instability. We consider the 
two-dimensional case and solve the Reduced MHD 
equations. [25j Denoting {0, ip} = d x (j)d y ip — d y (f)d x ip, 
we have 



d t vi0+{0,vi0} = {v,viv>}, 

d t ip + {(t>,il)} = rjVl^ + Eo. 



(4) 
(5) 



Here <fr and ip are the stream and flux functions of the in- 
plane velocity and magnetic field, so u = (~d y (p,d x (p), 
B = (—d y ip,d x ip). Eq. 0$ is the curl of the (inviscid) 
equation for an incompressible conducting fluid, Eq. ([S]) 



is the induction equation ([T]) uncurled. Eq is the equi- 
librium electric field, which must satisfy 8 x Eq = 0. The 
model of equilibrium flows described above corresponds 
to 4>q = Tgxy for \x\ < xq and (po = zLToXoy for \x\ > xq, 
which satisfies Eq. The magnetic-field profile ([3]) sat- 
isfies Eq. (J5]), provided we choose E = — va^q^g&ol for 
|ar| < x and Eq = — vaTqXq for |x| > Xq. 

Let us consider small perturbations to this equilib- 
rium, so ip = %jjQ + Sip, 4> = <Po + S<f>, an d linearize 
Eqs. (HHS]). If we seek solutions in the form S(p(x, y, t) = 
4>i{x, t) exp[ik(t)y] and Sip(x,y,t) = ipi(x, t) exp[ik(t)y], 
where k(t) = fco exp(— Tot), then (pi and ip\ satisfy 



(d 2 -k 2 )dt(Pi-T xd x (d 2 -k 
= [B Qy (x)(d 2 x - k 2 ) - 

d t ipi - Tox^Vi - B 0v (x)ik(j) 1 



)fa + 2T k 2 ( p 1 
B'^x^ik^, (6) 
= V (d 2 x -k 2 )^.(7) 



We now seek exponentially growing solutions, <pi{x, t) = 
— i$(x) exp(-ft) and ip\{x,t) = ^(x) exp^t). A solvable 
eigenvalue problem for 7 can be obtained if we assume 
7 ^> Tq , so the terms proportional to To can be neglected 
and kit) w fco. Note that the presence of the linear in- 
and outflows in this approximation is only felt via the 
equilibrium profile Bo y {x). Rewriting Eqs. ([5H7]) in terms 
of the dimcnsionlcss variable £ = x/5cs and denoting 
k = k v A /T = fcoics/2, e = (vT ) 1/2 /v A = 25 CS /L CS , 
and A = j/Tok, we get 

A($" - K 2 e 2 $) = -/(£)(*" - K 2 e 2 *) + /"(£)* ,(8) 
A*-/(e)$ = -(*"-K 2 e 2 *), (9) 

K 

where the derivatives are with respect to ^. 

The eigenvalue problem given by Eqs. jHHS]) is 
mathematically similar to the standard tearing mode 
problem, 22j except the role of resistivity is played by 
1/k. We shall assume that this parameter is small, i.e., 
we shall look for high-wavenumber perturbations of the 
current sheet. Another small parameter is e = (2/S) 1 ^ 2 , 
which is the inverse aspect ratio of the sheet. We shall 
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assume that kc <C 1 ■ All these assumptions will prove to 
be correct for the fastest growing modes. 

We now proceed as in the standard tearing mode cal- 
culation, considering first the outer region, £ ~ 1, and 
then the inner region, £ -C 1. 

Outer region. The behavior here is "ideal." Assuming 
that (1/k) < A < 1, we get from Eq. © * = A*//(f) 
and then from Eq. ([5]) 



no 



/(£) 



n 2 e 2 



(10) 



Given the functional form of the equilibrium /(£), solving 
this equation exactly is difficult. Instead, we shall solve 
it perturbatively, using K 2 e 2 <C 1. Neglecting « 2 e 2 to 
lowest order, we have an equation whose one solution is 
/(£). It is then easy to find the second solution, so the 
general solution can be written as follows 



dz 



(ii) 



where and arc constants of integration and ± 
refers to the solution at positive and negative values of 
£. We ask for the solution (but not its derivative) to be 
continuous at £ = 0. At small £, we have /(£) ~ a£, so 
the integral in Eq. (jlip is dominated by its upper limit 
and we find that Ct = = -a*(0). 

The constants C 1 are found by matching the solution 
(jTTJ) to the outer solution at |£| > £o (outside the current 
sheet). There we have, instead of Eq. (TIT)]) , 



= K 2 e 2 *, 



(12) 



so K 2 e 2 can no longer be neglected. The solution that de- 
cays at £ — * ±oo is *i! ± = exp(q=Ke£). Matching this 
solution and its derivative to the solution and using 
/(±f ) = ±1 and /'(±£o) = 0, we get C± = ±a*(0)/ K e 
and = a*(0)exp(Ke^ o )/K£- Thus, for |£| < £ 0) 

^^ 0^(0) ^ ^ & 

K e 7± ?0 / 2 (£) 



and for |£| > £ , 



(14) 



This outer solution is plotted in Fig. [3l It has a discon- 
tinuous derivative at £ = and is a classic unstable eigen- 
function known from the tearing-mode problem. The in- 
stability parameter A' = [*'(+0) - *'(-0)]/*(0) is 



KC 



So 



(15) 



where [. . . ] denotes the nonsingular part of the integral. 
The second term in Eq. (TIB")) , which is equal to 0.57, can 
be dropped compared to the first, so to lowest order we 




FIG. 3: The outer solution evaluated for e 1 = 10 3 at k 
340, corresponding to the maximum growth rate. 



have simply A' = 2a 2 /Ke.|3ll] 

Inner region. Here ( <C 1 and we can assume /(£) 
a£. Eqs. [9]), assuming 3> 1, become 



A3>" = -a£tf", 



(16) 
(17) 



Since A' = 2a 2 / ne 1, this eigenvalue problem is math- 
ematically similar to the one solved by Coppi et qZ. [26| 
for large-A' tearing modes. It reduces to solving the fol- 
lowing transcendental equation for the growth rate 



r ((A 3 / 2 - l)/4) 
r((A 3 / 2 + 5)/4) 



A'-^. (18) 

KC 



where T is the gamma function and A = Aa~ 2 / 3 K 1//3 . The 
width of the inner region is given by S = (A/k) 1 / 4 ^ -1 / 2 . 
Recall that A = 7 /T k and k = IcqVa/Fq- 

Eq. (|18p has two interesting limits: assuming A <C 1, 



(19) 



7 /r «1.63 K - 2 / 5 e- 4 / 5 ; 
assuming A — * 1— (from below), 



7 /r » ( aK ) 2 / 3 



3a 



(20) 



The maximum growth rate 7 max lies between these two 
asymptotics. Comparing the two terms in Eq. (|20[) shows 
that it is attained for ft max ~ e -3 / 4 ^> l and, therefore, 
7max/ro ~ e -1 ^ 2 3> 1. We note also that K max e ~ e 1 / 4 , 
A ~ e 1 / 4 and the inner layer width 8 ~ e 1 / 4 for the fastest 
growing mode. This confirms all of the ordering assump- 
tions we made in our calculation. 

Fig. [4] shows the dependence 7 (k) resulting from the 
numerical solution of Eq. (fT8|) for two different values of 
the current sheet aspect ratio e _1 . Both scalings (TT9"|) and 
([20]) are manifest. The vertical lines identify the values 
of k for which ne = 1. Strictly speaking, our calculation 
is only valid for values of k significantly to the left of this 
line. This is, however, of secondary importance, as the 
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FIG. 4: (Color online) The growth rate calculated from 
Eq. (|18|1 for two values of the aspect ratio: e _1 = 10 3 and 
10 6 . 

maximum of the growth rate lies well within the region 
of validity of our asymptotic ordering for e _1 > 10 3 . 

We have shown analytically that SP current sheets 
with large aspect ratios are intrinsically unstable to high- 
wave-number perturbations. Since e _1 = (5 f /2) 1 / 2 , the 
maximum growth rate scales with the Lundquist num- 
ber as 7 m ax ~ S 1 ' 4 ^. Typical Lundquist numbers in 
plasmas of interest are extremely large (e.g., S w 10 12 in 
the solar corona) , so the instability of the current sheet 
is extremely fast compared to the Alfven time ~ . 
One immediate consequence of this is that stable current 
sheets with aspect ratios above some critical value can- 
not exist. Numerical simulations [H|,[l6| suggest that this 
value is ~ 10 2 , corresponding to S c ~ 10 4 . Above this 
value, the sheet breaks up and a chain of plasmoids is 
formed. Their number scales as K max ~ S 3 ' s and their 
length is ~ 5" 1 / 8 larger than the current-sheet width <5cs- 

We note that, as far as we know, no numerical simu- 
lations of current sheets with aspect ratios significantly 
exceeding ~ 10 2 have been reported. At these val- 
ues, our asymptotic theory is not yet rigorously appli- 
cable, so we cannot test it against the existing simula- 
tions of SP reconnection. A nonasymptotic extrapola- 
tion of our results suggests that only a very small num- 



ber of plasmoids is to be expected — these are, indeed, 
seen. [H M EE EE H3, HI 

Since the growth rate we have found is much larger 
than the inverse Alfven time va/ '-^cSj the plasmoid width 
should become comparable to the inner-layer width 5 ~ 
S~ 1 / S 5cs before the plasmoids can be expelled from the 
current sheet by the Alfvenic outflows. The plasmoid 
evolution then becomes nonlinear. In order to have a 
quantitative theory of how the plasmoids affect the re- 
connection, one needs to understand what happens in the 
nonlinear regime. The dynamics in this regime will be 
dictated by a competition between three processes: the 
nonlinear growth and saturation of the plasmoids due to 
reconnection, the plasmoid coalescence, [H, [29[ and the 
expulsion of the plasmoids along the current sheet by the 
Alfvenic outflows. Numerical results on the stalling of 
the coalescence instability at large S suggest that multi- 
ple plasmoids can survive in the nonlinear regime rather 
than coalescing into a single plasmoid. [2 7j The formation 
of multiple large plasmoids has also been reported. [l8l[30j 
If the width of the plasmoids can indeed become bigger 
than the width of the current sheet before they coalesce 
or are expelled, estimates of the SP reconnection time 
must no longer be based on the parameters of the origi- 
nal current sheet but rather on some effective reconnec- 
tion region whose width is determined by the saturated 
plasmoid chain. A nonlinear study of the current sheet 
instability and plasmoid dynamics, as well as comparing 
the theory proposed in this Letter against more sophisti- 
cated and richer physics models, can only be performed 
by means of numerical simulations. We are presently de- 
veloping a computer code tailored to the study of this 
problem. Results will be presented in a separate publi- 
cation. 
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